
*********************************************
*** Program to Produce Figure 3 of appendix *
*********************************************

clear
drop _all
matrix drop _all
set more off

cd "$US_Ineq_Repl"
local Processed = "$US_Ineq_Repl/Processed"
local FFL = "$US_Ineq_Repl/Scripts/FFL"

// customize your file path
cd "`FFL'"
global data "`c(pwd)'\Data"
global figure "`c(pwd)'\Figures"

*** import data and variables ***

global expvar union pubsect manuf nonwhite female partte married smsa ///
	exper exper2 expe2 expe3 expe4 ///
    educ2 educ3 educ4 educ5 edex ee2 ee3 ee4 ee5 ee6 ee7 ee9 ee10 ee11 ee13 ee14 ee15 ///
    indu2 indu3 indu4 indu5 indu6 indu7 indu8 indu9 indu10 indu11 indu12 indu13 indu14 ///
	indu15 indu16 indu17 indu18 indu19 indu20 ///
	state2 state3 state4 state5 state6 state7 state8 state9 ///
	state10 state11 state12 state13 state14 state15 state16 ///
	state17  state18 state19 state20 state21 state22 state23 ///
	state24 state25 state26 state27 state28 state29 state30 ///
	state31 state32 state33 state34 state35 state36 state37 ///
	state38 state39 state40 state42 state43 state44 state45 ///
	state46 state47 state48 state49 state50 state51

use $data\us86.dta , clear
append using $data\us15.dta
append using $data\us86rw15.dta
keep lrwage3 time $expvar orgwgt	// use imputed wages if top-coded
compress
order lrwage3 $expvar
save $data\temp012.dta , replace

do runite.do

*Create Files to store iterations
forvalues t=5(5)95 {
	global qt=`t'
	local f="$data\"+"$qt"+"tempfig3.dta"
	
	capture confirm file "`f'"
	
	if c(rc) {
		di "File `f' does not exist. Running First iteration.."
		capture matrix drop Rt4qtau
		runite 1 $qt
		
		clear
		svmat Rt4qtau, name(b)
		rename b1 Dy
		rename b2 Dx
		rename b3 Dw
		rename b4 EXcov
		rename b5 EXoth
		rename b6 EXeduc
		rename b7 EXocc
		rename b8 EXind
		rename b9 EXres
		rename b10 UNcov
		rename b11 UNoth
		rename b12 UNeduc
		rename b13 UNocc
		rename b14 UNind
		rename b15 UNcons
		rename b16 UNres
		rename b17 TOTcov
		rename b18 TOToth
		rename b19 TOTed
		rename b20 TOTocc
		rename b21 TOTind
		rename b22 qtau
		rename b23 rep
		gen TOTDx=Dx+EXres
		gen TOTDw=UNcov+UNeduc+UNocc+UNoth+UNind
		replace qtau=round(qtau, 0.01)
		gen qt=round(100*qtau) 
		tostring qt, replace force
		sort qtau rep

		compress
		save "`f'" , replace
	}
}

local bsps=1000

forvalues t=5(5)95 {
	global qt=`t'
	local f="$data\"+"$qt"+"tempfig3.dta"
	
	use "`f'", replace
	quietly: sum rep, d
	local max_rep=`r(max)'
	if `max_rep'<`bsps' {
		local lb=`max_rep'+1
		forvalues j=`lb'(1)`bsps'{
			global ite=`j'
			do oneite.do 
			append using "`f'"
			compress
			save "`f'" , replace
		}
	}
}

*append all simulations
local qt=5
local f="$data\"+"`qt'"+"tempfig3.dta"
use "`f'", replace

forvalues t=10(5)95 {
	local f="$data\"+"`t'"+"tempfig3.dta"
	append using "`f'"
}

sort qtau rep

*create PEs and CIs
capture matrix drop ME
capture matrix drop MEt
capture matrix drop MEqt
capture matrix drop SEqt
capture matrix drop PE
capture matrix drop LL
capture matrix drop UL

forvalues qt=5(5)95 {
	di "doing quantile `qt'"
	
	foreach var of varlist D* EX* UN* TOT* {
		di "`var'"
		preserve
		keep if qt=="`qt'"
		quietly: sum `var'
		matrix ME=(nullmat(ME),1.96*sqrt(`r(N)'/(`r(N)'-1))*`r(sd)')
		restore
	}
	
	matrix MEt = ME,`qt'
	matrix drop ME
	matrix MEqt=(nullmat(MEqt)\MEt)
}

matlist MEqt

keep if rep==1
local vls = ""
foreach var of varlist D* EX* UN* TOT* {
	local vls= "`vls'" + "`var' "
}

local vls= "`vls'" + "qtau"
di "`vls'"
keep `vls'
order `vls'
gen qt=round(100*qtau) 
drop qtau

quietly:ds
mkmat `r(varlist)', matrix(PE) 

matrix LL=PE-MEqt
quietly:ds
local llvs = ""
foreach var of varlist `r(varlist)'{
	local llvs= "`llvs'" + "L`var' "
}
matrix colnames LL=`llvs'

matrix UL=PE+MEqt
quietly:ds
local ulvs = ""
foreach var of varlist `r(varlist)'{
	local ulvs= "`ulvs'" + "U`var' "
}
matrix colnames UL=`ulvs'

matrix ToT=PE,LL,UL
clear
svmat ToT, names(col)
drop Lqt Uqt
gen qtau=qt/100


*** plot ***
cd "$figure/coef"
set scheme s1color
global option1 xlab(0(0.2)1.0) ylab(-0.1(0.05)0.25) yscale(range(-0.12 0.2)) ///
	yline(0.0, lw(thin) lc(black) lp(dash)) ///
	xti("Quantile") yti("Log Wage Change", margin(0 2 0 0))
global option2 pos(11) ring(0) col(1) region(lstyle(none)) symxsize(8) ///
	keygap(2) textwidth(25)

twoway (connected Dy qtau, m(d) clw(medium) mc(red) lc(red)) ///
	(rarea LDy UDy qtau,col(red%10)) ///
	(connected Dw qtau, m(t) lp(dash) clw(medium) mc(green) lc(green)) ///
	(rarea LDw UDw qtau,col(green%10)) ///
	(connected Dx qtau, m(o) lp(shortdash) clw(medium) mc(blue) lc(blue)) ///
	(rarea LDx UDx qtau,col(blue%10)) ///
	, $option1 subtitle("A. Change in Log Wages 2015-1986") ///
	legend(order(1 3 5) lab(1 "Total Change") lab(3 "Structure") lab(5 "Composition") ///
	$option2 ) saving(fig3a,replace)


graph twoway (connected Dx qtau, m(o) lp(shortdash) clw(medium) mc(blue) lc(blue)) ///
	(rarea LDx UDx qtau,col(blue%10)) ///
    (connected TOTDx qtau, m(d) lp(dash) clw(medium) mc(teal) lc(teal)) ///
	(rarea LTOTDx UTOTDx qtau,col(teal%10)) ///
    (connected EXres qtau, m(t) clw(medium) mc(dkorange) lc(dkorange)) ///
	(rarea LEXres UEXres qtau,col(dkorange%10)) ///
    , $option1 subtitle("B. Aggregate Composition Effects") ///
    legend(order(1 3 5) lab(1 "Composition") lab(3 "Explained") lab(5 "Spec. Error") $option2) ///
    saving(fig3b,replace)


*** combine and export ***
graph combine fig3a.gph fig3b.gph , xsize(7) ysize(10) col(1) saving($figure\fig3.gph, replace)
graph export $figure\fig3.pdf, replace

graph export $US_Ineq_Repl\Results\Figures\FigF3.png , replace

rm fig3a.gph
rm fig3b.gph
rm $data\temp012.dta
